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iv.  A  computational  method  has  been  developed  in  which  adjacent  cells  in  a  cell  netwoik  interact  through  an  exchange  of 
paiticles.  commensurate  with  the  equilibriiun  fluxes  of  mass,  momentum,  and  eneigy.  This  corresponds  to  the  inflnite 
collision  rate  limit  of  the  model  gas,  resulting  in  very  low  viscosity.  The  scheme  is  more  efficient  than  any  direct  simulation 
with  cellular  autmnata  or  other  method  for  solving  the  lattice  Boltzmann  equations. 

V.  A  simple  multi-speed  model,  the  nine-velocity  model  has  been  studied  extensively.  Solving  the  (inviscid)  Reimann  problem 
with  the  model  yields  almost  all  of  the  phenomenology  associated  widi  a  perfect  gas.  An  exaa  (viscous)  shock  profile  is 
computed  for  the  model  and  is  compart  to  a  Navier-Stokes  shock  profile.  Adiabatic  channel  flow  is  simulated  with  the 
model  and  the  results  compared  to  an  integral  solfltion  of  the  htavier-Stokes  equation.  The  comparisons  in  both  cases  are 
excellent.  It  is  also  shown  that  the  nine-velocity  gas  does  not  permit  steady  supersonic  flow. 
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Investigations  of  the  Motion  of  Discrete- Velocity  Gases 


1.  Statement  of  Work 

1.  Develop  the  IDSMC  for  engineering  use.  including: 

i.  Introduce  a  more  realistic  molecular  model. 

ii.  Implement  diffuse-wall  boundary  conditions. 

2.  Investigate  othermethodologies  for  simplifying  the  treatment  of  molecular  collisions  in  the  IDSMC. 

3.  Bench  mark  the  IDSMC  with  a  large-scale  computation  of  engineering  importance,  for  example, 
perfect-gas  flow  over  a  sphere  or  cylinder. 

4.  Continue  investigations  of  relaxation  to  equilibrium  and  anisotropy  in  the  2-D  9-velocity  cellular 
automaton.  Consider  the  applicability  of  dynamical  systems  theory  to  questions  of  reversibility  and 
stability. 

5.  Extend  the  2-D  9-velocity  CA  to: 

i.  3-D  15-vclocity  model. 

ii.  2-D  25-velocity  model. 

Hi.  3-D  45-velocity  model. 

6.  Extend  the  2-D  CA  to  three  dimensions. 

7.  Compute  test  cases  with  the  above  models  for  comparison  with  the  IDSMC  and  DSMC. 

2.  Status  of  the  Research  Effort 

The  accomplishenmts  of  this  grant  ate  described  in  the  publications  listed  in  Sec.  6  and,  more  recently,  in 
the  Caltech  Ph.D.  Thesis  “A  Study  of  Multi-Speed  Discrete  Velocity  Gases”  by  Balu  Nadiga.  They  include: 

1.  Development  of  a  direct  simulation  Monte  Carlo  method  for  rarefied  gasdynamics,  patterned  after 
the  DSMC  of  G.  A.  Bird,  in  which  molecular  velocities  ate  discretized  to  be  integers  (IDSMQ.  The 
method  has  beat  implemented  on  coarse-grained  multi-computers  such  as  the  Intel  iPSC,  Symult 
2010  and  Intel  Gatiuna  with  a  self-adtqrtive  rectangular  mesh  to  optimize  load  balancing.  Many 
problems  in  gasdynamics  were  solved  as  test  beds  and  to  measure  performance. 

2.  Formulation  and  implemention  of  a  flnite-dilferrace  method  for  solving  the  discrete-velocity  (lat¬ 
tice  gas)  Boltzmatm  equations.  Calculations  were  made  of  problems  in  heat  transfer,  shock  wave 
structure  and  vs^r  dqrosition,  typically  with  9-  or  13'Velocity  models,  but  in  one  case  with  a  40- 
velocity  model.  As  with  all  of  our  numerical  direct  simulations  of  discrete-velocity  gases,  most  of 
the  repetitive  calculatimis  ate  dmie  tmly  cmce,  and  the  results  are  stored  in  look  up  tables. 

3.  Study  of  the  iq^licid>ility  of  multi-speed  discrete-velocity  gases  to  compressible  flows  from  a  fun¬ 
damental  point  of  view.  The  equation  of  state,  the  anisottr^ies  and  the  advection  velocities  for 
multi-speed  models  on  the  square  and  triangular  lattices  are  derived.  Dqrendence  on  the  model  of 
any  of  these  to  leading  orrter  in  the  flow  velodiy  is  shown  to  be  only  through  a  fourth  rmnnent  of 
the  stationary  equilibrium  qreed  distribution.  Thus  the  pathologies  shown  by  multi-speed  lattice 
^tses  dm  be  made  rather  small. 
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4.  Development  of  a  computational  schone  (EFM)  in  which  adjacent  cells  in  a  cell  netwotk  interact 
through  an  exchange  of  particles,  commotsurate  with  the  equilibrium  fluxes  of  mass,  momentum, 
and  oieisy.  This  corresponds  to  the  infinite  collision  rate  limit  of  the  model  gas,  resulting  in  very 
low  viscosities.  The  scheme  is  more  efficient  than  any  direct  simulation  with  cellular  automata  and 
other  methods  for  solving  the  lattice  Boltzmann  equations. 

5.  Exhaustive  study  of  a  simple  multi-speed  model,  the  nine-velocity  model.  Solving  the  shock  tube 
problem  with  the  model  yields  almost  all  phenomenology  associated  with  a  perfect  gas.  An  exact 
shock  profile  is  computed  for  the  model  and  is  compared  to  a  Navier-Stokes  shock  profile.  Adiabatic 
channel  flow  is  simulated  with  the  model  and  the  results  compared  to  an  integral  solution  of  the 
Navier-Stokes  equation.  The  comparisons  in  both  the  cases  are  excellent.  It  is  also  shown  that  the 
nine-velocity  gas  does  not  permit  steady  supersonic  flow. 

In  the  following  sections  we  describe  results  which  have  not  yet  been  published. 

3.  A  Mettiod  for  Near-Equilibrium  Discrete- Velocity  Gas  Flows 

We  have  developed  a  simulation  schone  for  discrete-velocity  gases  based  on  local  thermodynamic  equi¬ 
librium.  Exploitation  of  the  kinetic  nature  of  discrete-velocity  gases  results  in  a  natural  splitting  of  fluxes, 
and  the  resultant  scheme  strongly  resembles  the  original  processes.  Calculation  of  exact  molecular  interac- 
timis  by  direct  simulation  is  computationally  intensive,  perhaps  uimecessarily  so,  and  limits  the  maximum 
achievable  Reynold’s  number.  It  is  felt  that  a  better  simuliaion  technique  for  discrete-velocity  gases  can  be 
achieved  by  using  in  conjunction  with  the  underlying  concept  of  particles,  an  Eulerian  picture  of  the  flow 
field.  The  resulting  scoiario  is  one  in  which  cells  in  a  flow  field,  constituted  a  discrete-velocity  gas,  inter¬ 
act  in  terms  of  mass,  momoitum,  and  energy  fluxes  at  the  boundaries,  the  currency  of  the  interaction  being 
the  particles  with  the  discrete-velocities.  Further,  if  the  interactions  between  the  cells  are  based  on  the  local 
equilibrium  fluxes,  then  we  have  a  method  which  is  kineflc  based  and  at  the  same  time  models  the  infinite 
collision  rate  (near-equilibrium)  limit  of  the  kinetic  ptxx^.  These  ideas  come  from  the  equilibrium  flux 
method  (EFM)  of  Pullin,  which  is  a  kinetic-theory-based  finite  volume  method  for  calculating  the  flow  of  a 
compressible  ideal  gas.  In  the  EFM,  the  particle  ensonble  simulation  of  an  ideal  gas  is  replaced  by  a  scheme 
in  which  adjacent  cells  in  a  cell  netwotk  interact  directly  through  exchange  of  mass,  momentum,  and  energy, 
the  exchange  in  t  being  calculated  using  the  equilibrium  flux  relations. 

The  kinetic  nature  of  the  scheme  and  the  modeling  of  the  irfinite  collision  rate  limit  result  in  a  low 
value  of  (numerical)  viscosity,  whose  macroscopic  behavior  is  qualitatively  like  that  of  a  viscosity  resulting 
from  molecular  interactions.  A  first  order  method  and  two  second  order  methods  using  the  total  variation 
diminising  principle  are  tteveltqred.  An  example  application  is  presented.  While  this  technique  retains  the 
attractive  features  of  lattice  gases  —  simplicity  and  parallel  evolution,  it  represents  a  large  improvement 
over  the  currently  used  lattice  gas  automata  simulations  and  lattice  Boltzmarm  equation  solution  techniques, 
in  that  the  evolution  is  based  on  the  inviscid,  non-heat  conducting  limit  With  this  ^proach  much  higher 
Reynolds  numbers  are  adiievatrie  than  with  any  present  methods,  given  the  same  cmnputer  resources.  The 
ideas  being  general,  the  scheme  is  applicable  to  any  discrete-velocity  model,  and  to  lattice  gases  as  well. 

4.  An  Exact  Shodk  Sointion  for  a  Nine- Velocity  Gas 

The  dK)dr-atructure  proUem,  because  of  the  absence  of  solid  boumlaries  and  the  simplicity  of  the  ge¬ 
ometry,  is  particuitfly  usdiii  in  studying  the  ifliysics  of  computatkoal  models  of  fluids.  The  success  of 
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disciete-velocity  models  of  fluids,  inclusive  of  a  class  of  lattice  gases,  in  r^ioducing  various  fluid  phenom¬ 
ena  is  now  well  known.  Most  woik  in  this  area,  however,  has  been  confined  to  single-speed  models,  where 
the  absence  of  m  independent  eneigy  variable  leaves  the  thermodynamics  incomplete.  While  the  exact  struc¬ 
ture  of  shock  waves  in  such  single-speed  models  are  known,  their  counterparts  in  the  nuUtiple-speed  models 
are  likely  to  be  more  important  because  of  the  non-trivial  thermodynamics. 

An  exact  shock  solution  is  computed  in  a  multiple-speed  discrete-velocity  gas,  the  nine-velocity  gas, 
the  simplest  multiple-speed  model  on  the  square  lattice  in  two  dimoisions.  The  multiplicity  of  speeds  en¬ 
sures  nontrivial  thermodynamics.  The  model  Boltzmann  equations  are  solved  numerically  for  steady  one¬ 
dimensional  flow  in  an  infinite  domain  by  tracing  out  ^proximate  solution  curves  in  three-dimensional  phase 
space  using  the  dux)ting  method.  This  is  one  of  the  few  instances  where  an  exact  soluticm  has  been  obtained 
for  the  nine-velocity  gas.  G>mparison  to  the  shock  structure  in  a  monatomic  perfect  gas,  as  given  by  the 
Navier-Stokes  equation,  shows  excellent  agreement  The  shock  in  the  nine-velocity  gas,  has  an  overshoot 
in  entropy  like  in  a  monatomic  gas.  The  near-equilibrium  flow  technique  for  discrete-velocity  gases  (see 
above),  a  kinetic  flux-splitting  method  based  on  the  local  thennodynamic  equilibrium,  is  also  seen  to  capture 
the  shock  structure  at  the  level  of  the  finite-difference  discretization  (cell  size)  remarkably  well. 

5.  2-D  Adiabatic  Channel  Flow  by  Direct  Simulation  with  a  Lattice  Gas 

The  viscous  and  compressible  but  laminar  flow  of  a  perfect  gas  through  an  adiabatic  charmel  is  not  a 
fully  understood  problem.  Its  incompressible  counterpart  constitutes  Hagen-Poiseulle  flow  and  is  an  exact 
solution  of  the  Navier-Stdces  equations.  There  have  beat  few  instances  of  a  full  simulation  of  the  problem 
of  laminar  compressible  viscous  charmel  flow.  One  of  the  problems  is  the  imposition  of  the  downstream 
boundary  conditions.  In  this  section,  a  simulation  of  the  channel  flow  using  the  nine-velocity  model  is 
described.  The  lattice  gas  method  being  a  particle  method,  the  imposition  of  the  boundary  conditions  is 
simpler. 

The  computational  domain  is  rectangular.  The  interaction  of  the  system  with  the  external  surroundings 
is  only  through  the  perimeter  of  the  domain.  With  flow  along  the  longer  dimoision,  the  interactions  along 
the  diorter  sides  constitute  the  inflow-outflow  boundary  conditions  and  the  interactions  along  the  longer 
dimoisions  describe  vriiat  htqrpens  on  the  charmel  walls. 

The  simulation  of  the  wall  boundary  condition  is  simple.  Two  possible  wall  collisions  are:  i)  specular 
reflection;  only  the  component  of  velocity  perpendicular  to  the  wall  is  reversed.  In  this  case,  the  accom¬ 
modation  coefficient  of  the  wall  is  zero,  i.e.,  flow  properties  parallel  to  the  wall  are  unaffected  by  the  wall. 
The  flow  velocity  at  the  wall  can  be  non-zero,  ii)  “bounce-back";  the  velocity  of  the  particle  impinging  the 
wall  is  reversed.  The  accommodation  coefficient  in  this  case  is  two.  This  is  a  rather  severe  case  of  a  wall 
from  which  the  particles  do  not  bounce  in  random  directions,  but  instead  ranonber  their  incoming  velocity 
and  reverse  it  on  collision.  A  realistic  wall  has  an  accommodation  coefficioit  of  unity  (diffuse  wall)  and 
this  is  easily  achieved  by  uniformly  but  randomly  choosing  between  the  two  types  of  reflection.  In  both 
the  ^Kcular  and  botmce-badc  reflection  of  particles  at  the  wall,  the  particle  speeds  are  preserved,  «isuring 
adidraticity  of  the  wall. 

The  inflow  boundary  condition  is  implemented  through  a  "reservoir,”  a  body  of  fluid  udtich  is  unaffected 
the  flow  in  the  diannel.  The  reservoir  in  the  present  case  is  a  small,  indqrendent  and  steady  state  con¬ 
stantly  updated  region,  periodic  in  the  flow  direction,  located  upstream  of  the  actual  dianrrel.  A  particular 
lattice  otrfunm  of  the  reservoir  is  assumed  to  be  located  immeditdely  upstream  of  the  diannel  entrance,  so 


that  panicles  at  this  lattice  column  directed  downstream  are  introduced  at  the  channel  entrance  during  the 
following  time  step.  In  the  computaticms  presented  here,  the  upstream  reservoir  uses  a  mobius-strip  bound¬ 
ary  OHidition  for  the  walls  parallel  to  the  flow,  and  a  periodic  boundary  condition  in  the  flow  direction,  as 
introduced  1^  Kadanoff  et  al.  in  1989.  The  flow  in  the  upstream  reservoir  is  sustained  by  imposing  a  body 
force,  the  level  of  which  is  adjusted  to  obtain  the  required  flow  rate.  This  generates  a  velocity  profile  which 
is  parabolic  for  low  forcing  levels. 

The  out  flow  boundary  condition  at  the  downstream  «ul  of  the  computational  channel  amounts  to  the 
impositimi  of  a  downstream  pressure  at  the  last  column  of  lattice  sites.  That  colunm  is  modeled  as  a  porous 
wall  with  the  porosity  depending  inversely  on  the  downstream  pressure  to  be  simulated.  For  a  given  pressure, 
the  porosity  of  the  wall  at  any  time  is  uniform  across  the  wall,  but  changes  location  at  each  time  step  to 
preclude  the  formation  of  any  upstream  patterns. 

The  multi-computer  is  configured  as  a  linear  array  of  processors  along  the  flow  direction.  The  compu¬ 
tational  domain  is  equally  partitioned  along  the  flow  direction  and  the  partitions  assigned  to  the  nodes  (one 
cmnputer  of  the  multi-computer).  There  are  no  regions  of  overly  among  the  various  nodes  and  the  intem- 
ode  communication  is  limited  to  moving  particles  either  leaving  or  entering  the  domain  of  a  node..  All  the 
cmnputations  presented  here  were  done  with  the  Intel  iPSC/860  GAMMA  parallel  machine,  which  consists 
of  64  nodes  of  iPSC/860  wired  in  a  hypercube  architecture.  The  cormnunications  in  the  present  case  were 
handled  under  the  Reactive  Kernel  environment  of  Seitz  et  al. .  The  sustained  lattice  site  update  rate  for 
these  computations,  with  the  present  implementation,  was  about  one  million  sites  per  second. 

Noteworthy  features  of  the  results  are: 

•  The  flow  accelerates  slowly  at  first,  but  much  more  rrqridly  at  higher  velocities.  The  last  S0%  of  the 
increase  in  velocity  is  over  about  12  widths  of  the  charmel  (the  full  channel  is  128  widths  long). 

•  The  velocity  profile,  initially  parabolic,  flattens  near  the  downstream  end. 

•  The  pressure  profile  across  the  channel  is  quite  flat  up  to  relatively  high  flow  velocities. 

•  The  flow  is  nearly  isothermal  (width-averaged  tonperature  almost  constant  in  the  downstream  di¬ 
rection)  at  low  velocities,  rather  than  constant  density,  as  might  have  been  expected. 

•  The  width-averaged  pressure  drops  off  nearly  linearly  down  the  channel  at  low  velocities,  becoming 
increasingly  non-linear  at  higher  velocities. 

•  The  width-averaged  density  and  pressure  profiles  are  very  similar  at  low  velocities. 

•  The  effects  of  viscous  dissipation  are  confined  to  a  region  near  the  walls.  A  rise  in  entropy  observed 
in  the  central  core  is  mainly  due  to  the  falling  density. 
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